

%% load estimation results
load('parameter_estimates_round2_replication.mat')
[~,use] = min(l0);
thetaStar = theta(:,use);
lambdaStar = theta((size(thetaStar,1)-5):end,use);
UStar =  payoffs(thetaStar,X,Z,YY,reb,mpow);
clear K l0 srngSV SV theta theta0 time use Y exitn;
theta = thetaStar;


N=size(X,1);
KZ=size(Z,2)/5;
% NS=size(V_R,2);
nBeta=length(reb{1});
nPhi=length(mpow{1});
nGamma=length(reb{2});
nChi=length(mpow{2});
nYY=size(YY,1);

Xrebel = X(:,reb{1}); Xmpow = X(:,mpow{1});

ZUS = Z(:,1:KZ);
ZUK = Z(:,KZ+(1:KZ)); 
ZFRN = Z(:,2*KZ+(1:KZ));
ZRUS = Z(:,3*KZ+(1:KZ));
ZCHN = Z(:,4*KZ+(1:KZ));

ZrebelUS = ZUS(:,reb{2}); ZmpowUS = ZUS(:,mpow{2});
ZrebelUK = ZUK(:,reb{2}); ZmpowUK = ZUK(:,mpow{2});
ZrebelFRN = ZFRN(:,reb{2}); ZmpowFRN = ZFRN(:,mpow{2});
ZrebelRUS =ZRUS(:,reb{2}); ZmpowRUS = ZRUS(:,mpow{2});
ZrebelCHN = ZCHN(:,reb{2}); ZmpowCHN = ZCHN(:,mpow{2});


gamma_US=theta(nBeta+4+nPhi+1,1);
gamma_UK=theta(nBeta+4+nPhi+2,1);
gamma_France=theta(nBeta+4+nPhi+3,1);
gamma_Russia=theta(nBeta+4+nPhi+4,1);
gamma_China=theta(nBeta+4+nPhi+5,1);
gamma0=theta(nBeta+4+nPhi+5+(1:nGamma),1);

disp(' ')
disp('mean of gammas:')
disp([mean(gamma_US + ZrebelUS * gamma0)...
mean(gamma_UK + ZrebelUK * gamma0)...
mean(gamma_France + ZrebelFRN * gamma0)...
mean(gamma_Russia + ZrebelRUS * gamma0)...
mean(gamma_China + ZrebelCHN * gamma0)])

disp(' ')
disp('std. dev. of gammas:')
disp([std(gamma_US + ZrebelUS * gamma0)...
std(gamma_UK + ZrebelUK * gamma0)...
std(gamma_France + ZrebelFRN * gamma0)...
std(gamma_Russia + ZrebelRUS * gamma0)...
std(gamma_China + ZrebelCHN * gamma0)])

